#include <iostream>
#include <fstream>
#include <functional>
#include <mpi.h>
#include "config.hh"
#include "utils/hexagon.hh"
#include "utils/solver.hh"
#include "utils/debug.hh"
#include "kernel/methods.hh"

void output(icehpc::hexagon const& hex, icehpc::solver const& s);

int main(int argc, char** argv)
{
    // 初始化 mpi
    MPI::Init(argc, argv);

    // 定义了一个 method map 来做 dispatch，目前只有 jacobi，所以下面也就直接写
    // 名字了
    std::map<
        std::string,
        std::function<void(icehpc::hexagon, icehpc::solver&)>
    > method{
        { "jacobi", icehpc::jacobi }
    };

    // 获取进程数 size 和当前进程的 rank
    auto size = MPI::COMM_WORLD.Get_size();
    auto rank = MPI::COMM_WORLD.Get_rank();

    // 解析配置文件
    auto toml = cpptoml::parse_file("config.toml");
    icesp::configuration::config config{toml};
    // division 就是报告中讲的将六边形水平的分成 4 * div 块，方便定点标号
    int    division{config.cases.begin()->second.division};
    // eps 控制迭代误差
    double eps     {config.cases.begin()->second.eps};

    // hexagon 管理正六边形相关的数据，以及一些工具函数
    icehpc::hexagon hex(division);
    // solver 管理每个 mpi 进程所用的数据
    icehpc::solver s(hex.n, size, rank);
    // 对 solver 进行初始化，主要是初始化其中的 s.b，stencil 中的常数项，也可以
    // 理解成一个线性方程组 Ax=b 其中的 b
    icehpc::init(hex, s);

    // 输出相关的信息，这个 PRINT3 在 debug.hh 中实现，可以通过定义
    // DEBUGGING_ENABLED 来使用
    PRINT3("division = ",     division,    " ");
    PRINT3("eps = ",          eps,         " ");
    PRINT3("active point = ", hex.nactive, " ");
    PRINT3("hex.n = ",        hex.n,       " ");
    PRINT3("hex.c = ",        hex.c,       " ");
    PRINT3("this is ",        rank,        "\n");

    // 主循环，不断进行迭代
    for (auto loop = 0; ; loop++) {
        // 调用 jacobi 进行计算
        method["jacobi"](hex, s);

        // 通过 mpi 的 Allgather 讲每个进程的局部解 s.u 同步到全局解 s.au 上
        MPI::COMM_WORLD.Allgather(&s.u.front(), s.tot, MPI::DOUBLE, &s.au.front(), s.tot, MPI::DOUBLE);

        // 通过 mpi 的 Allreduce，每个进程都计算出最大误差
        double error;
        MPI::COMM_WORLD.Allreduce(&s.max_error, &error, 1, MPI::DOUBLE, MPI::MAX);

        // 输出迭代次数和当前误差
        if (rank == 0) {
            PRINT3("iteration ", loop, " ");
            PRINT3("error ", error, "\n");
        }
        // 如果误差达到要求就退出
        if (error <= eps) break;
    }

    // 0 进程负责输出
    if (rank == 0) output(hex, s);

    MPI::COMM_WORLD.Barrier();
    MPI::Finalize();
}

void output(icehpc::hexagon const& hex, icehpc::solver const& s)
{
    std::ofstream fout{"output/points.dat"};
    // output to tecplot data format
    fout << "Title=\"solution\"\nVariables=\"x\", \"y\", \"u\"\n";
    for (auto i = 0; i < hex.n; i++) {
        if (hex.outside(i)) continue;
        auto p = hex.real_pos(i);
        fout << std::fixed << std::setprecision(6)
            << p.x << " " << p.y << " " << s.au[i] << "\n";
    }
}

